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Abstract 

>' 
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JiJ ■ The ferroelectric polarization of perovskite oxides is much larger than im- 

plied by displacement of static ionic charges. We use an explicitly correlated 
0^ . scheme to investigate the phenomenon; charge transport is evaluated as a 



geometric quantum phase. Both covalence and electron-electron interaction 



i-rt ■ enhance polarization in the weakly correlated regime. At higher values of the 

O ■ electron-electron interaction, the system undergoes a transition from a band 



insulator to a Mott insulator: the static ionic charge is continuous across the 
transition, whereas the polarization is discontinuous. Above the transition, 
oxygen transports a positive charge. 
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The large values of the spontaneous polarization in ferroelectric perovskites cannot be 
explained on the basis of a simple ionic model. For instance, the effective dynamical charge 
associated to an oxygen displacement during the ferroelectric distortion is of the order of 
—6, i.e. about three times the nominal static value given by a completely ionic picture 
|T|J^ (atomic units of charge are adopted throughout). First-principles calculations, per- 
formed within the local-density approximation (LDA), have demonstrated that this effect 
is due to a large amount of covalence |P] : as a basic feature of such calculations, the effect 
of electron-electron interaction is accounted for in a mean-field way. In the present paper 
we investigate a different facet of the problem, namely how the picture is changed by ex- 
plicitly considering correlated wavefunctions: such investigation is based on a completely 
different theoretical and computational scheme. In the moderately correlated regime, we 
find that the electron-electron interaction reduces the static ionic charges but — in apparent 
contradiction — significantly enhances the macroscopic polarization. At higher values of the 
electron-electron interaction the system has a transition from a band insulator to a Mott 
insulator, first investigated by Egami et al. Q]. When we examine the transition from the 
standpoint of macroscopic polarization, several exotic phenomena show up. In particular the 
static ionic charges are continuous, whereas the polarization is discontinuous: therefore the 
polarization (or equivalently the associated geometric quantum phase [^,^) is the primary 
order parameter for the transition. 

The ferroelectric perovskites have formula unit ABO3, where A is a mono- or divalent 
cation, and B is a penta- or tetravalent cation in the transition series. The simplest ferroelec- 
tric structure is the tetragonal one, where the most important symmetry-breaking distortion 
is a relative displacement of the B sublattice with respect to the O one. The essence of the 
phenomenon can be schematized by considering only the linear 0-B chain along the fer- 
roelectric axis, which is symmetric in the paraelectric (so-called prototype) structure, and 
asymmetric in the ferroelectric one. In order to make a correlated wavefunction available, we 
schematize the Hamiltonian of such a linear model ferroelectric with a two-band Hubbard 
model at half filling, first introduced by Egami et al. 0], and whose macroscopic polarization 



has not been investigated so far. In order to reduce the number of parameters, we arbitrarily 
set equal values of the Hubbard U on anion and cation. The Hamiltonian of the prototype 
structure is therefore 



El (-l)'^^^- - io(49+i, + H.c.) ] + Uj^nn-n^l , (1) 
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and depends on two parameters besides U: the hopping to, and the difference in site energies 
Eb — Eq = 2 A. We restore charge neutrality by assigning a classical charge of +2 to the 
cationic sites: this ensures that in the extreme ionic limit (to = and U < 2A) the total 
static charges are ±2. The energy difference Eb — Eq is of the order of 3 or 4 eV in the 
materials of interest 0, taking e.g. Ti or Nb as the B cation; it also equals the gap of the 
noninteracting model. In the following of this study we assume A = 2 eV, and to = 3.5 eV, 
such as to get a mean-field band width of 5.3 eV, which is a realistic value for the valence 
bands of these materials. The fact that to is of the order of 2A clearly indicates a mixed 
ionic/covalent character. 

We plot in Fig. [I|(a) the static ionic charge as a function of U (triangles). The calculations 
have been performed in 8-sites supercells, using skew (quasiperiodic) boundary conditions 
on each electronic variable, and then taking the average over the boundary conditions (alias 
over the supercell quasimomenta k), as in Ref. P]. This ensures that the U ^ limit is 
numerically equal to the fully converged thermodynamic limit of the noninteracting calcu- 
lation. Typically we have used 30 k points, corresponding to 120 k points in the unfolded 
Brillouin zone of the noninteracting system. The numerical diagonalization was performed 
via the well known Lanczos algorithm, which provided the ground state electronic wave- 
functions \E'o(/i;) with an energy tolerance close to the machine accuracy (lO^^^to ), with 
less than 100 Lanczos iterations in all cases. In order to reduce the problem size we have 
explicitly used the conservation of the number of spin-up N-^ and spin-down A^^^ particles. 
The subspace with Ni = N^ = 4, where the ground state lies, contains only 4900 elements, 
thus allowing a large reduction of the full Hilbert space (amounting to 4^ = 65536 elements). 
Use of the center-of-mass translation symmetry is also possible, but becomes useful only for 



larger systems. In an insulating system as the present one the size effects are small, and 
further minimized by fc-averaging ^ . For a few parameter values we have indeed performed 
12-sites calculations, and checked that the results agree with the 8-sites ones to within one 
per cent, at least for the quantities studied here. 

The ferroelectric distortion is the zone-center optical phonon, where the site coordinates 
are displaced by a relative amount ^, typically of the order of ^f = 0.05a, where a is the lattice 
constant: this distortion asymmetrically affects the hopping matrix elements. We assume 
the other parameters fixed during the ferroelectric distortion, while for the t variation we 
assume the simple Su-Schrieffer-Heeger [Q linear dependence t = to ± 2aC,- For reasons 
given below, a realistic modeling of these materials requires a rather large electron-phonon 
coupling a, typically aa ~ 10 eV in the mean-field case. The tight-binding noninteracting 
Hamiltonian is trivially diagonalized as 



e{k) = ± V A^ + Atl cos2 ka/2 + 16(02 sin^ ka/2 . (2) 

The band structure is quadratic in ^, hence the (linear) deformation potential vanishes in the 
prototype structure. Nonetheless the band shift induced by the actual ferroelectric distortion 
^F is rather large (about —0.8 eV at the zone boundary) in semiquantitative agreement with 
LDA calculations (see e.g. Fig. |l](a) in Ref. 0). 

The static charges are somewhat reduced by the distortion, as shown in Fig. |l|(a), circles. 
When U is increased to large values, the system undergoes an interesting transition, from a 
band insulator to a Mott insulator. With our parameter values, the transition occurs at Uc = 
2.27to- Egami et ai, who first investigated this transition using only k = wavef unctions, 
found a discontinuous drop in the static ionic charges of the prototype structure, while the 
charges of the ferroelectric structure were found continuous as a function of U. We reproduce 
their results, but we also find that that the discontinuity disappears as fc-point convergence 
is approached. A careful analysis is displayed in Fig. [l|(b), which incidentally proves the 
effectiveness of the k average |Q] in getting rid of spurious finite-size effects. The apparent 
discontinuity is due to a level crossing which occurs at A; = and not at k ^ 0, as discussed 



below. We have explicitly verified that the computed discontinuity is inversely proportional 
to the number of k points used. Furthermore the discontinuity disappears even with a coarse 
mesh if the mesh is displaced on the k axis in order to avoid the k = singular point: this 
is also shown in Fig. |l|(b), open circles. 

We are interested in the macroscopic polarization AP induced by the distortion, when 
the sites are continuously displaced from the prototype structure (^ = 0) to the ferroelectric 
one (^ = ^f)- We therefore need evaluating how much charge is transported along the 
chain during a relative displacement of the two sublattices, in a vanishing electric field P| ; 
if we choose to keep the origin fixed on a cationic site, the transport is purely electronic. 
The electronic charge transport is best evaluated as a geometric quantum phase, as first 
shown — for the explicitly correlated case — by Ortiz and Martin . The rationale behind the 
geometric phase approach is that the dynamical charge is a quantum-mechanical current — 
hence a phase of the wavefunction — and bears in general no relationship to the modulus of 
the wavefunction when periodic boundary conditions are used. 

The starting ingredients are again the same ground-state wavefunctions \E'o(/^) discussed 
above. One then removes the Bloch-like phase upon defining $o(^) = '^o{k)Y\exp{—ikxj), 
where Xj are the site coordinates. The wavefunctions $o(^) so obtained are periodic over 
the supercell at any k, and implicitly depend on the Hamiltonian parameters. Only the 
dependence upon ^ and U will be relevant in the following discusssion. One then defines the 



many-body generalization of the geometric phase first introduced by Zak [jlO 
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{^,U)=zJdk{^o{k)\^Mk)), (3) 



where the fc-integration is over the Brillouin zone of the supercell. The numerical calculation 
proceeds as in the uncorrelated case |^ , and the macroscopic polarization of the ferroelectric 
structure is 

AP=[7(eF,f/)-7(0,f/)]/27r, (4) 

defined modulo a polarization quantum of magnitude 1, corresponding to the transport of 



one charge over one cell. This is one half of the quantum of the mean-field theory |^, where 
double occupancy of one-particle states is enforced. 

In the special case where the electron-phonon coupling a is taken as vanishing, then each 
site may only transport its static ionic charge (shown in Fig. |I|): in fact the polarization 
calculated as a geometric phase accounts precisely for this rigid charge transport. Notice 
however that the two alternative calculations are not nnmenca//?/ equivalent, thus providing 
a useful convergence test. Using our typical fc-point mesh given above, the error is smaller 
than 10~^. When a 7^ 0, the dynamical charge is no longer equal to the static one, and is 
typically much larger than it. 

Let us first illustrate the mean-field (f/ = 0) calculation. We take aa = 10 eV, which 
provides a spontaneous polarization in agreement with both measurements and LDA calcu- 
lations in typical ferroelectric perovskites. The polarization AP is almost linear in ^: the rel- 
evant quantities to display are therefore the average dynamical charge (Z*(^)) = aAP{^)/^ 
and the linear (or Born) dynamical charge of the prototype structure Z* = aP'(O). The 
static (cationic) charge is 1.47: charge transport is enhanced by a factor larger than four by 
the electron-phonon coupling, thus providing giant dynamical charges, and large values of 
the spontaneous polarization. The actual values within our model are {Z*{C,f)) = 5.95 and 
Z* = 7.28. 

We now switch discussing the effects of electron-electron interaction. The first interesting 
phenomenon occurs already in the centrosymmetric system, having real wavefunctions ^o{k): 
the Zak phase 7 changes discontinuously by vr at the transition point. Equivalently, one finds 
a Berry phase of vr around the rectangular loop in the {k, U) plane shown in Fig. |^, since the 
vertical sides of the rectangle do not contribute 0J^. This means that the real wavefunction 
undergoes a sign change when transported along the closed path: the commonest occurrence 
of such a feature, well known in molecular physics |]11|, is due to the presence of a point 
of double degeneracy inside the domain. This is precisely the case here: there are two 
well distinct states whose energies cross at the point (0, Uc), whereas at fc 7^ there is no 
degeneracy. We have numerically checked the level crossings by exploiting the metastability 



of the Lanczos iteration across the transition. Coming now to the physical meaning of such 
transition, we notice that a phase change of n corresponds to the transport of an electronic 
charge over half a lattice constant, from an oxygen site to a cationic one. Notice once more 
that such transport occurs without affecting the static charges. This is a virtue of the ring 
geometry of our chain, whereas in a finite linear chain, owing to continuity, charge transport 
would obviously affect the static charges of the end sites. We thus discover that the geometric 
phase — and not the static ionic charge — is the primary order parameter for the transition. 



One could even straightforwardly generalize the "band-center operator" of Ref. ||I0[ to the 
many-body case, and characterize the transition by saying that the crystalline ground state 
is an eigenstate of such operator belonging to different eigenvalues below and above Uq- 

We then consider the polarization of the ferroelectric structure, calculated as in Eq. (^ 
and where the $o(^) ^^re complex. We plot in Fig. ^ the average dynamical charge (Z*) as 
a function of U, for several values of ^. In the moderate-?/ region below Uc the electron- 
electron interaction produces a significant enhancement of the polarization. Notice that this 
latter feature is in apparent contradiction with the fact that the static ionic charges decrease 
instead with increasing U. The most prominent feature of Fig. |^ is the divergence at Uq, 
which has an interesting physical meaning. The (near) divergent curve, corresponding to our 
smallest ^, is an approximation to the Born dynamical charge: the figure then indicates that 
at U = Uc an infinitesimal sublattice displacement, (starting from the prototype structure) 
induces a finite charge transport, hence an infinite Z*. At finite ^ values instead — and in 
particular at the value C,f = 0.05a corresponding to a realistic ferroelectric distortion — the 
polarization has a finite and large discontinuity at Uc- Notice that the Zak phase of the 
ferroelectric structure is continuous as a function of U, and therefore both the divergence 
and the discontinuity of the dynamical charges must be traced back to the discontinuity of 
the centrosymmetric [C, = 0) Zak phase in Eq. (^). 

In conclusion, we have investigated here the effect of electron-electron interaction in 
the macroscopic polarization of ferroelectric perovskites, by means of an explicit model 
Hamiltonian which captures the basic features of the phenomenon, and exploiting the geo- 



metric phase approach. The very large polarization of these materials owes to their mixed 
ionic/covalent character, and is further enhanced by electron-electron interaction, as long as 
its strength remains moderate. At higher strength, the system undergoes a transition from 
a band insulator to a Mott insulator. At the transition point, the polarization is discontinu- 
ous, and even reverses its sign for a given sublattice displacement. In the highly correlated 
regime the cation transports a negative dynamical charge, and oxygen a positive one. 
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FIGURES 
FIG. 1. Static charge of the cation as a function of the Hubbard U. (a) Prototype centrosym- 

metric structure (triangles) and ferroelectric structure (circles) below the transition, (b) Different 

computations for the prototype structure at the transition point Uc- Empty triangles: 10 k points. 

Filled triangles: 100 k points. Empty circles: 30 k points, displaced such as to avoid A; = 0. 

FIG. 2. Average dynamical charge of the cation as a function of the Hubbard U, for different 
values of the displacement ^. In order of increasing value of the discontinuity at Uq, the curves 
represent: C = Cf = 0.05a; ^ = 0.035a; ^ = 0.0245a; ^ = 0.014a; ^ = 0.0035a. 

FIG. 3. Rectangular loop in the (k, U) plane which encircles the level crossing at /c=0 and 
U=Uc, in arbitrary units. The projection over the k axis coincides with the Brillouin zone of the 
supercell. The cross is at the degeneracy point. 
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